function [ skel ] = mat2( img, tolerance )

%Binarization
bw = blockBinarization(img, 1);
[M, N] = size(bw);
bw(1,:) = 1;
bw(M,:) = 1;
bw(:,1) = 1;
bw(:,N) = 1;
bw = imfilter(bw, fspecial('gaussian', 5));

% Distance Transform with euclidean distance
dist = bwdist(bw);

% Filtering with a LOG
tmp = imfilter(dist, fspecial('average', 5));
tmp = imfilter(tmp, -fspecial('log'));
figure('name', 'log'); imshow(tmp, []);

% Normalization
tmp = tmp - min(tmp(:));
tmp = tmp / max(tmp(:));
figure('name', 'normalization'); imshow(tmp);

% Binarization with optimal threshold
tmp = im2bw(tmp, tmp(1,1) * (1 + tolerance / 100.0)) .* (1-bw);
figure('name', 'threshold'); imshow(tmp, []);

% This morphological operator removes isolated pixels
% tmp = bwmorph(tmp,'majority', 1);

eroded = imerode(tmp, ones(5));
tmp = imreconstruct(double(eroded), double(tmp));
figure('name', 'reconstruct'); imshow(tmp, []);

% Thin the skeleton for get 1 pixel width
thin = bwmorph(tmp,'thin', Inf);
dist = double(dist) / max(dist(:));
figure('name', 'result'); imshow(bw * 0.8 + dist .* 0.2 - tmp .* 0.2 + thin);

skel = thin;

end